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We study the spin-pseudospin Hamiltonian of the Ising Model in Transverse Field (IMTF) for 
^ ■ pseudospins, coupled to the XY-spins on a triangular lattice. This model appears from analyses 

' of the quarter-filled ladder compound NaV20s, and pseudospins represent its charge degrees of 

freedom. In the molecular-field approximation we find that the model possesses two phases: charge- 
^ . disordered without spin gap; and a low-temperature phase containing both the anti-ferroelectric 

(zigzag) charge order and spin dimerization (spin gap). The phase transition is of the second 
kind, and the calculated physical quantities are as those one expects from the Landau theory. 
One of particular features of the phase diagram is that the inter-ladder spin-pseudospin coupling, 
responsible for the spin gap generation, also destroys the IMTF quantum critical point, resulting in 
the exponential behavior of T c in the region of Ising's coupling where the IMTF is always disordered. 
We conclude that our mean-field results give a qualitatively correct description of the phase transition 
in NaV205, while a more sophisticated analysis is warranted in order to take into account the thermal 
fluctuations and, probably, the proximity of the IMTF quantum critical point. 
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I. INTRODUCTION 

Properties of the models with coupled spin and orbital degrees of freedom have been actively studied from early 
70-th, mostly in the context of the Jahn- Teller transition metal compounds^ In recent years it has been a growing 
theoretical effortSiiiSiSiLSi^iiS in order to understand the ground-state properties and excitation spectra of a particular 
(high-symmetric) class of one-dimensional spin-orbital models. The latter can also be viewed (up to some slight 
modifications at most) as two spin- jj Heisenberg chains coupled via a bi-quadratic exchange. This interest is motivated 
by various unusual states, occurring in many kinds of low-dimensional magnetic systems: e.g., TDAE-Cgo^ two- 
leg spin-i ladder compounds (for a review, see Ref.fTlj): recently discovered spin-gapped materials Na2Ti 2 Pn 2 
t^- [ (Pn = As, Sb)i£ and NaV^Os? 1 ^ which do not fit into the standard spin-Peierls scenario. 

For the one-dimensional (ID) spin-orbital model there is a special point (i.e., coupling's values) where the effective 
Hamiltonian has an extra SU (4) symmetry, and the exact solution has been known for quite a long time^* Another 
exactly solvable special point for the SU(2) (8 5(7(2)-symmetric Hamiltonian was found recently by Kolezhuk and 
Mikeskai^ Apart of these special values of couplings, the available results on the lD-spin-orbital models are known from 
mean-field analyses, renormalization group, bosonization, conformal field theory, numerics, and/or their combinations. 

As was re-emphasized recently^ from the view point of application to real materials we are mostly interested in, 
high symmetry of a spin-orbital Hamiltonian is rather artificial. Even if the total spin is conserved from physical 
■ grounds, i.e., the symmetry of the Hamiltonian's spin sector is indeed SU(2), the symmetry of the orbital sector is 
usually lower (e.g., Ising or XY-type). So, applicability of the results obtained from high-symmetric Hamiltonians 
(e.g., SU{2) £g> SU (2)) to the cases of low-symmetric Hamiltonians following from microscopic considerations of a given 
Q ' problem, is not so straightforward^ 

, It has been argued that the spin-gapped compounds Na2Ti2Pn2 4 and NaV^Ogii can be modeled by a two-band 
quarter-filled Hubbard Hamiltonian. In the regime of strong on-site Coulomb repulsion such Hamiltonian can be 
mapped onto a spin-orbital modelji so the appearance of a spin gap is due to interplay between the charge and orbital 
degrees of freedom. Instead, one can work directly with the electronic Hamiltonian. Phase diagram of the ID quarter- 
filled two-band Hubbard model (Hubbard ladder) has been recently studied both numerically^ and analytically^ 

In this paper we will explore a particular scenario for the spin gap mechanism in a spin-orbit model, where the orbital 
( "pseudospin" ) degree of freedom corresponds physically to the ordering charge displacement (disproportionation). 
Such a model appears in the applications to NaV^Osi^S The quarter-filled ladder-system 1 ^ NaV^Os is up-to-date the 
only well-established transition-metal compound with a localized spin- i moment distributed, in the high-temperature 
phase, equally over two transition-metal ions in a mixed valence state. This unique situation of two I/ 4 - 5 ions sharing 
one spin s = i leaves room for the charge disproportionation 2V 4 ' 5 — > V 4 +V 5 occurring below the critical temperature 
of T c = 34 K, together with the opening of a spin-gap (without a long-range order for spins). At T c the unit cell 
doubles both along the ladders (6-direction) and along the rungs (a-direction), see Ref.[2lJ for a review. 
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There have been several mechanisms proposed for the driving force of this phase transition: the spin-phonon 
coupling like in a spin-Peierls system^ the electron-phonon coupling, the intra-ladder Coulomb repulsion between V 
and I^ 5 -ions on the neighboring rungs^SLSiiS^, or a combination of those interaction terms^ The experiments seem 
however to rule ouiM^Sk the original proposal^ for the spin-Peierls scenario, while the zigzag charge ordering pattern 
predicted by the Coulomb-driven scenarios of this phase transition)22iS4i2^ is now established crytallographically. 30 
Another argument in favor of the key role played by the charge ordering at this phase transition, comes from the fact 
that the coupling between spins and zigzag-ordered charge leads to2£ the observed dispersion^ along the a-direction 
(perpendicular to the ladders) of the gapped magnon mode in the low-temperature phase. 

Since NaV2 0s is an insulator, one can model the charge degrees of freedom in NaX^Os on a given rung by a 
pseudospin T. The system of quarter-filled coupled ladders can be mapped onto a triangular lattice where one spin 
and one pseudospin reside on the same site22i25i2&, compare Fig.^ Mostovoy and Khomskii in their study of NaV205 
proposed the following mechanism for the spin dimcrization driven by charge ordering (i.e., spin gap generation in 
the charge ordered phase) the antiferroelectric (AFE) in-ladder ordering (i.e., zigzag charge ordering) is due to 
the inter-rung Coulomb repulsion. The ordering pattern is shifted by a half of the inter-rung distance from ladder 
to ladder (see FigQ, and it results in alternation of the spin-exchange coupling along the ladder induced by its left 
and right neighbors. Then, as a consequence of the spin-exchange dimerization, the systems acquires a spin gap. 
Mostovoy and Khomskii^ also pointed out that for this effect to occur, it suffices for the dimerization term of the 
spin-exchange coupling (i.e, the inter-ladder spin-pseudospin coupling) to be linear over the charge order parameter. 
This mechanism resembles so the familiar spin-Peierls scenario, with the proviso that phonons are replaced by the 
charge order parameter. 

Going deeper into this suggestion for the coupling which could be a key to the problem of NaV2 0s, we perform a 
detailed analysis of such a phase transition starting from a simplest possible Hamiltonian compatible with symmetry 
restrictions. The effective Hamiltonian is that of the Ising Model in Transverse Field (IMTF) for pseudospins T, 
coupled to the spin-i Ay-chains. Since the phase transition in NaV20s is thermal in nature, we treat the charge 
(pseudospin) sector of the Hamiltonian in the molecular-field approximation, while retaining in our analysis the 
experimentally observed ID nature of the spin excitations. Recent x-ray«2ii2£ and NMRS experiments showed large 
fluctuation regions on the both sides of T c in NaV20s, and the order parameter critical index is found to be ~ 
0.17 — 0.19. It is close to (3 = 1/8 of the 2D Ising model, indicating along with the correlation lengths measurements^ 
on the 2D-Ising universality class of this phase transition. According to Ref.JHI] the ID character of structural 
fluctuations manifests itself only at T > 60K. This allows us to expect the mean- field treatment to be a qualitatively 
reasonable approximation for the problem we are dealing with. 

We find that our spin-pseudospin model has a phase diagram with two phases: the disordered one without charge 
order or spin gap; and the ordered low-temperature phase with both the anti-ferroelectric (zigzag) charge order 
and spin gap. One of our most important findings is that the model always orders (continuously, via a second- 
order transition) into that phase at low temperatures. The critical temperature of this phase transition shows various 
regimes of dependence on model's parameters, but it never vanishes as far as the inter-ladder spin-pseudospin coupling- 
is non-zero. Physically, it means that the role of the inter-ladder coupling is twofold: this coupling not only creates 
the spin gap when the phase transition is driven by the charge (pseudospin) ordering in the IMTF, but it also destroys 
the quantum critical point of the IMTF and generates the charge-ordered and spin-gapped phase in the case when 
charge ordering in the pure IMTF is impossible. 

The rest of the paper is organized as follows. In Section^we define the model Hamiltonian and describe the methods 
employed to handle it. In Section IlIII we present the mean-field equations which govern the critical behavior of the 
system. Section ITVl presents the analysis of the phase diagram, i.e., the critical temperature of the anti-ferroelectric 
phase transition as a function of the Hamiltonian's couplings. The results for the order parameter, spin and pseudospin 
susceptibilities, specific heat capacity calculated in the vicinity of the anti-ferroelectric phase transition are given in 
Section The summary of our results and the discussion on their application for NaV205 are presented in the final 
Section EH 
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II. MODEL AND FORMALISM 



A. Model 



We choose the following effective Hamiltonian: 

H= — ^^Xnn + ^ ^mn%n,n+l + ^ S mn S m ,„+i (f) 

m.n m,n m,n 

+ SmnSm,n+l [^ST^„ n ^m ; ti+l + ^(^m+l,n+l ~ ^m-l,n)] 

where the spin operators for a given site are defined as S i = \a l (i — x,y,z) via the Pauli matrices a % . The 
pseudospin operators are related in the same way to the Pauli matrices t 1 : T l = \t % ■ We use the distinct notations 
for the Pauli matrices a, r in order to emphasize the fact that they operate in different spaces (spin, pseudospin). So 
the components, e.g., a x and t x do not have the same quantization axis. By labelling the two sets of Pauli's matrices 
cr, r with x, y, z we mean only their conventional representation and commutation relations. The subscript 1 < m < M 
stands for the number of ladder, while 1 < n < M numbers the rung. We assume the Hamiltonian coupling constants 
to be positive. Note in making comparison with the microscopic Hamiltonians derived for NaV20sj22i2£i2& that in our 
Hamiltonian we made the ^-rotation in the pseudospin space T x i— > T z , T z >— > — T x . 

The first two terms on the right hand side of Eq.Q correspond to the Ising Model in Transverse Field (IMTF), the 
third term to the Heisenberg Hamiltonian, and in the last one which describes the spin-pseudospin interaction oc e, 
we retain only linear term over the rung charge displacement operator T x . This inter- ladder coupling term is allowed 
by the point symmetry group D 2 h (the NaV20s space group is — Pmmn)^ Note that a similar intra-ladder 
spin-pseudospin coupling J2 mn ^mn^m,n+i(T^ n+1 — T£ n ) would be odd with respect to the mirror-plane through 
the bridging-oxygen (the 6c-plane), therefore it is forbidden by symmetry. 

In the following we will work with the dimensionless Hamiltonian Ti. = H/Q 

T~£ = ~ y ] %nn "t" 7^9 ^""^ %nn%n,n+l (2) 
m , n m , n 

+ Dmn [J + ^%nn%n : n+1 + £ {%n+l,n+l — %n~l,n)] I 

m,n 

dimensionless temperature T — > T/f2, and dimensionless couplings 

g = Ji/Q, J = J H /^, (3) 
A = j ST /n , e = e/n , 

where the dimerization operator D mn = S mn • S m ,«+i- For the applications we have in mind, the range of the model's 
parameters under consideration will be restricted to the cases 

(J> A) < 9, e < (g, max{ J, A}) (4) 
In the present study we consider the Hamiltonian with the XY spin-spin interaction, i.e., 

n — n XY — q x q x i qy qv rr.\ 

J-^mn ^mn — u mn u m,n+l ' u mn u m,n+l \ u ) 

= — ( Q+ Q~ _1_ <?~ Q+ 

2 ^ mn u m,7i+l u mn m,n+l 

where we used the conventional definition S l± = S x ± iS y . 



B. Molecular Field Approximation 



We proceed with the approximate treatment of our model in the following way: The model Hamiltonian © can 
be viewed as H — Himtf + Hxy- "^imtf is the Hamiltonian of the IMTF in terms of the pseudospin operators T 
which is exactly solvable by its own. Dependence of the effective coupling J^ n on T in the XY spin Hamiltonian 
TtxY = X) JmnDmn precludes the total Hamiltonian from being exactly solvable. 
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We apply a version of the Mean Field Approximation (MFA) which is rather known as the Molecular Field method 
(Approximation) in the theory of phase transitions 3 ^. This approximation assumes separability of the total Hamilto- 
nian density matrix p over S and T, i.e., 

p = p s ®p T (6) 
and the following single particle (operator) ansatz for the pseudospin density matrix 

P T = Il^n. (?) 

pLn = -^-exp(-/3h m „T m „) (8) 

mn 



Here h m „ is the on-site vector of the molecular field to be defined self-consistently, (3 is inverse temperature (we set 
the Boltzmann constant ks = 1), and 

^ n = 2cosh*^ (9) 

in the on-site partition function. 

Since the independent spin-pseudospin averaging © with the single-particle density matrix JSJ) makes J e ff a function 
of the average values of T, no further approximation is needed for the density matrix p s . The statistical mechanical 
problem with 

p s = ^-e- f3Us , Z s = Tre~' 3Ws , (10) 

Ti-s = Ti-xy (here Tr includes only spin-operator states) is exactly solvable^ and will be dealt with in the next 
subsection. To summarize more qualitatively, we solve exactly the problem of independent spin XY chains with their 
exchange constants determined by the molecular fields of quasi- two-dimensional pseudospins. 



C. Spin Hamiltonian 

The standard Jordan- Wigner transformation 3 ^^ (JWT) when applied independently to one-dimensional spin oper- 
ators on a given ladder, results in fermionic operators which commute on different ladders (see, e.g., Ref. |37j). Strictly 
speaking this does not cause problems for our choice of Hamiltonian J2J j since its spin part does not couple different 
ladders. However we find this situation with the ID JWT somehow unpleasant, especially keeping in mind possible 
generalizations of the present Hamiltonian for future work. There are generalizations of the JWT for d > 1 available 
in the literature (see, e.g., Ref. 38 39]). 

We will employ a version of the 2D JWT proposed by Azzoua 3 ^, particularly convenient for our model. We also 
include to the spin Hamiltonian a term generated by the uniform external magnetic field 

H-cxt = ^ hMS mn (11) 

The transformation is given by the following equations^ 

S mn = exp(-i7r ^ n ki)c mn (12a) 

fc,2G2l m „ 

Stm = 4mexp(i7T n kl) ( 12b ) 

where n mn = <} mn c mn is the fermionic number operator, and the above summations include all rungs and ladders 
lying on the left from the m-th ladder plus rungs lying below n-th rung at the m-th ladder, i.e., 

k, I e %nn <=>{l<k<m-l,l<l<Af}U{k = m,l<l<n-l} (13) 

The fermionic operators satisfy the canonical anticommutation relations 

{cL.*} = dmkSni (14a) 

{Cmn,C k l} = (14b) 
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and 

D mn = 2 ( C lnn c m,n+1 + c m,n+l c "in) (15) 
S-mn = c lnn c mn ~ (^) 

For the effective coupiing J^ n determined by the mean-field values of the pseudospin operators we take a dimerized 
ansatz, so 

^XY — ^ i^mn ) Jmn = A m + ( — 1)"-Bm (17) 

Applying then the JWT defined above and a Bogolyubov transformation^ we obtain 

Hs = Wxy + W^t = - X -MMh M + ^ E mv {q)d) mqv d mqv (18) 

where the new fermionic operators d mqu also satisfy the canonical anticommutation relations, and their spectrum is 



E„ w (q) = h M + v\j A 2 n cos 2 qa + sin qa (19) 

The extra index v = ±1 is due to dimerization, i.e., doubling of fermion species, summation over q includes the 
reduced Brillouin zon e [— and a is the distance between rungs (we will set a = 1 in the following). For more 
details, sec, e.g., Rcfs. 40 41] . The operator of the total z-component of spin per rung reads 

m,n m,g,i/ 

The spin contribution to the free energy per rung from the Hamiltonian 7ig (|18f) [cf. Eqs.^UJ] is 

In 2 \ r% (5E mu {q) 



(3 -K/3M 



f dq In cosh ; 



III. MEAN-FIELD EQUATIONS 

The molecular fields are defined from the condition that they minimize the mean-field free energy T = (Ti.) + T (In p) . 
The angular brackets stand for averaging with the mean- field density matrix © , and the Hamiltonian Tt is given by 
Eq.©. We consider the case when the external magnetic field Km is absent. (Note that contrary to /s l|21ll . the 
quantity fj- defined in the same fashion via Z T (j5J is not the pseudospin free energy). 

We take the following ansatze for the Ising pseudospin magnetizations (i.e., the charge ordering parameters in terms 
of the real physical quantities) 

= m z (22) 

(TZn) = (-l) m+n m x (23) 

So, similar to the case of the IMTF (g > 0) we assume the possibility of the anti-ferroelectric (AFE) in-ladder charge 
ordering (i.e., the zigzag ordering^) , alternating however from ladder to ladder. It is easy to see from the Hamiltonian 
(|2*|) that ansatz (|23[) creates a dimerization in the spin sector, so a natural assumption for the dimerization operator 
average is 

(D mn )=-[t+{-l) m+n 8\ (24) 
Then the average mean-field energy per rung in the thermodynamic limit M.N — > oo is 

£mf = = -m z - ^gm 2 x - (Ji + \tm 2 z + 2eSm x ) (25) 
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The molecular Weiss fields h z mn = -h z ; h x mn = (-l) m+n+1 h x ; J^ n = a + (-l) m +™& are defined by the following 
equations: 

h z = -^F = f + 2Xtm z (26a) 
om z 

h x = - ^ MF — gm x -f 2e5 (26b) 
om x 

a = _^I = J + Am 2 (26c) 

b = d&MF = 

do 

The order parameters I|22I23I24[1 obtained as partial derivatives of the corresponding quantities /s,/r with respect 
to their conjugate Weiss fields iJ^SJ, should be determined from the system of four coupled equations 

f 1 + 2\tm z J3t> 

m z = tanh — 27a 

m x = 9 + ^ £f] tann (27b) 

t = - ( 2 d^^-^tanh^(v) = -tn(Aj) (27c) 

= — / 2 dV^-r tanh /3£M s —rtniAj) (27d) 

where f) = y/h x + fv* is the absolute value of the Ising molecular field, and instead of 5 we use the new parameter rj 

5 = m x r] (28) 

We also introduced the auxiliary parameters 



£(ip) = y cos 2 + A 2 sin 2 <^9 (29a) 

2em x , „, , 

A = 29b 

J + Am 2 

= f(^ + Am 2 ) (29c) 



The following useful relationship holds at any temperature 

m 2 +m 2 = itanh 2 ^ (30) 
Alternatively, the mean-field equations (|27|l can be derived from the free energy (per rung, Iim = 0) 

/ = 2\tm\ + ^gm 2 x + 2e5m x - Tin 4 (31) 
rflncosh^ 1 2 



— / dip In cosh f3£ (ip)) , 

7T In ' 



understood as f(m z ,m x ,t,S), by minimization over its variables, [f) is defined according to Eas. (|26|l .] Thus, / plays 
a role of the Landau functional, defining parameters m z , m x , t, 8. 

Before presenting the analysis of Eas. (|27|) . let us describe their main properties more qualitatively. The system 
of these four coupled equations determines four unknown mean-field parameters m z ,m x ,t,r) (5) as functions of the 
effective Hamiltonian couplings and temperature. The first couple of Eqs.J37|) in the case A = e = is familiar from 
the mean-field treatment of the pure IMTF (see, e.g, Ref . |34j ) . and has some similar properties with the pure case 
even for non-zero A and e. There is always a non-trivial solution < m z < \ for the field-induced magnetization 
(the lower bound value m z — is attained in the limit g — > oo). Below certain critical temperature T c a non-trivial 
solution m x ^ appears. m z varies continuously across T c , contrary to the case of the IMTF (A = e = 0) when it 
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stays constant (= 1/g) for all T < T c . Non-zero m x results in generation of the alternating term b in the exchange 
interaction J^n an d the non-zero dimerization parameter 6. Consequence of b ^ is a spin gap. Note that in the 
disordered phase (m x = 0), Ea. (|27t is an empty statement at T ^ 0, since 5 — as well. The parameter t is 
non-critical: it is some continuous non-negative function. 

At T < T c one equation from the pair H27b .b) can be written in the form 



m z = g + 



4e 2 



7r(J + \m 2 z ) 



—t T < T 

IT 



(32) 



more convenient for further analyses. 



IV. CRITICAL TEMPERATURE 



At T = T c we have Eqs.(|27K) as 



m z = - tanh — (H t n ) , 



(33) 



another equation for m z l)32[l. and parameters t n ,rj„ are given by Eqs.(j2Zt,d) with A = 0. The latter two functions 
have the following expansions: 



t„(0,x) 



x(l - \x 2 ) + C(x 5 ), x< 1 



4~V- 4 

1 _ 2_J_ 

L 24 x 2 ' 



X > 1 



and 



n (0 x)~ 1^ X{1 ~ ^ + 0(x5) ' X<1 



where A = 7re f_ 7 w 1.6685, and 7 » 0.5772 is Euler's constant. 



(34) 



(35) 



Case £ = 0: re-entrance 



In order to better understand the model, let us start with the analysis of the mean-field equations for the case when 
the dimerization coupling is absent, i.e, e = 0. Then from Eas. H32l33|) we find the critical temperature T c = 1//3 C 



m z g 



In 



l+2m z 
l-2m, 



where m z is the solution of the following equation 



_i 2A / J + Arn^ l + 2m 2 

m z =9 *n 0, -T. ln 1 5 — 

7r \ 2m z <7 1 — zm 2 

It can be shown from Eas. (|32l33l34|l that for large Ising coupling g ^> g\ where 
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2 1 



(36) 



(37) 



(38) 



the critical temperature T c sa g/A is large , while m z w 1/g is small. So, the behavior of these parameters is the 
same as in the pure IMTFi^ It can also be shown that at g = g\ the critical temperature T c becomes zero (m z = -|). 
Thus g\ plays a role of the renormalized (due to the coupling A) mean- field critical point g* = g\(X = 0) = 2 of 
the pure IMTF, where T c vanishes and no ordering in m x occurs for g < g*~ (The mean-field g* is two times less 
than the value known for the exactly solvable ID IMTF. For 2D and 3D IMTF the values of g* are known only from 
approximations and/or numerics, but they are bigger then the mean-field predictions^) 

However, the role of A is more subtle than a simple shift in the parameters comparatively to the pure IMTF. The 
critical line T c (g) manifests re-entrant behavior in the region g < g\ or, to state it differently, the function g(T c ) is 
non-monotonic, contrary to the case A = 0. (See Fig. 0). Also, T c (g\) — does not mean that g\ is the minimal 
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coupling g m in at which -ordering can occur. How pronounced the re-entrance is, i.e., the width of the re-entrance 
region 

5min < g < gx, (39) 

depends on the relative values of model's couplings J, A, g. A more detailed analysis of the re-entrance in this model 
will be presented in a separate paper 4^ At this time, we say that any value of the coupling A 7^ always generates the 
re-entrance on the phase diagram g(T c ), i.e., g m in — > g\ —> 2 only when A — > 0. Characteristic results of the numerical 
solution of our equations are given in Fig. At J <C A this feature is pronounced most [the case J = 0, A = 1 shown 
in this figure], while at J ~ A it is more smeared [see the curve at J = A = 1]. Notice that when g — ► g\, the curve 
T c (g) approaches zero normally to the g-axis. 

More qualitatively, the origin of the re-entrance can be understood from the free energy (|31f) which is minimized 
by our mean-filed equations. First two terms on the r.h.s. of Ea. l|31[) augment the free energy / with the increase 
of m z and/or m x , while at the same time the last two terms on the r.h.s. of that equation, explicitly proportional 
to the temperature, decrease / via the parameters t),/3, A. The interplay of this contributions to / involves, apart of 
the temperature, the couplings g, J, A. The latter possess a range [the re-entrance region <|39|l ] wherein the minimal 
free energy is achieved by re-distribution of the values of m z (T) and m x (T). In particular, within this region, upon 
say, decreasing T, the system goes smoothly from the disordered phase m x (T) — to the ordered one m x {T) 7^ 0, 
and back to m x (T) = 0. Note also from Ea. (|33|) . that A creates a self-consistent renormalization of the external field 
[fi, in terms of the dimcnsionful Hamiltonian l[T|l]. resulting in, particularly, the temperature dependence of m z in the 
ordered phase. (Compare to the pure IMTF where at T < T c it is frozen at m z = l/g). According to Ea. (|34|l the 
effective external field 1 + 1\t n m z ji: grows at low temperature, bringing m z closer to its maximal value 1/2 and in 
the same time, forcing m x to decrease [cf. Ea. <|r}()jl ] . This is another argument in the way to qualitatively explain why 
the system can re-enter the disordered phase upon decreasing the temperature. However, we should point out that 
these results and arguments are based on the mean-field analysis, so whether or not the re-entrance survives a more 
sophisticated treatment of the model, is not clear at this time. 



B. Case e / 0: double re-entrance 



An important qualitative change in the critical behavior of the model occurs upon switching on the dimcrization 
coupling e. At any e 7^ the critical temperature is never zero, even if g <C g\. Let us note first that when {A, T} — > 
then r) n ~ — ln(min{A, T}) is divergent, so this limit should be treated with care. When g 3> g\ [cf. also conditions 
the behavior of the critical temperature and magnetization m z is the same as in the case e = 0. Upon decreasing g, 
the critical temperature T c decreases, and according to Eq.(j23 m z « |, while from Eas. (|32l34l35() mj 1 « c\ — c%-\sxT c , 
where c\, C2 are some constants. From the consistency of those equations we conclude that for any finite e the critical 
temperature never vanishes. To put it differently, coupling e destroys the quantum critical point of the IMTF. This 
constitutes a very important feedback effect of the spin chains on the charge degrees of freedom. 

Characteristic numerical results for T c (g) are shown in Fig. El Similar to the case e = 0, the critical temperature 
shows the re-entrance in the region g < g\. Moreover, instead of going to zero at g = g\ as for e = 0, the curve T c (g) 
turns left, and T c remains finite even at g = 0, albeit exponentially small in e. For the two curves discussed above in 
the case e = 0, we show in Fig.|2]their counterparts at e = 0.1. As one can see from the whole curve T c (g), at bigger 
values of A the system manifests a well-pronounced double re-entrance. 

It follows from Eas l|32l33l34l35|) that the left turn of the critical temperature and its failure to vanish at g = gx 
can be analytically described as a BCS-type solution for T c (g), generated by finite e. Approximate solutions for T c (g) 
found in two regions of g are 



where 



T c « I . _ (40) 

I ^exp[ -g(g x - g)] , BCS regime 

J^J+\ (41) 

The boundary where the low-temperature BCS regime sets in and the related formulas are applicable, is given by the 
condition 

4e 2 

BCS regime : g < g\^ = (42) 

nJ 
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The field-induced magnetization m z in these regions is 



, 9 > 9X 

, BCS regime 



(43) 



C. Case e / 0, A = 0: no re-entrance 

The analytical treatment in the intermediate regime when g ~ g\ is involved due to re-entrance. The situation 
simplifies when A = and the re-entrance is absent (see the curve for J = 1, A = 0, £ = 0.1 in Fig. [21 ■ Then both the 
intermediate g > g* (g* = 2) and the BCS regimes are well described by the approximate equation 



* 4e 2 , AJ 
Q — Q H in — — 



4e" 



-l/Tc 



(44) 



In particular, at the IMTF quantum critical point on finds that it is destroyed by the inter-ladder dimerization 
coupling £, resulting in 



T c a In" 1 



ttJ 



> .9 = 9' 



(45) 



So when the re-entrance is absent, upon decreasing g the critical temperature T c monotonically decreases from the 
IMTF linear regime T c cx g via the logarithmic dependence (|4B|> towards the exponential BCS regime (14011 . It is 
interesting to note that a similar inverse-log dependence of the transition temperature near quantum criticality has 
been found in a recent molecular-field study of Cu2Te20sBr2, a spin-system containing coupled spin-tetrahedra™. 



V. PROPERTIES OF THE AFE PHASE 



A. AFE order parameter 



In the ordered phase we can derive from Eas. H26l27l30(l the following equation 



2 yj m| + ml \ 2 



: tanh 



y/ ml + ml ) = 1 



where 



9e = 9 + 



4e 2 



tt(J + Xm 2 z ) 



'In 



(46) 



(47) 



To determine the temperature behavior of the AFE order parameter m x in the immediate vicinity of the critical 
temperature where 



T c -T 

ST* 1 ' 



(48) 



we need to expand Ea. (|46|l near T c , taking into account that all parameters m x ,m z ,g £ entering this equation are 
temperature-dependent and related via the mean- field equations. 
It can be shown that to leading order 



t n (Aj) « t n (0j) + A 2 i nA 2 
r] n {Aj) m ?7„(0,/3) + A 2 r] nA 2 

with the partial derivatives given by the following equations 

dt n (Aj) 



^n, A 2 
Vn,A 2 



OA 2 

dVn(Aj) 



dA 2 



A=0 



A=0 



a(P)-t n (0,(3) 



l a 0) + l b 0)-r) n (OJ) 



(49) 
(50) 



(51) 
(52) 
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and 



<P) = 9 / — ; (53a) 

£ Jo cosh pz z 



HP) 



13 dz 



o cosh 2 (3z 



VI- z 



2 



2/3 



-tanh/3z 



(53b) 



Using the above equations in the leading-order expansion of Eas. (|32l4fi() near T c , we can show that 

m x cc t (54) 

along the whole line of the critical temperature depicted in Fig. [21 i.e., T c defines a second-order phase transition, and 
the order parameter has the mean- field critical index 1/2. In general, the coefficient of proportionality in the above 
relationship has to be defined numerically. For the regimes of strong Ising coupling and BCS, the order parameter 
can be calculated analytically. 

For the functions defined by Eas. l|53|) the following expansions are obtained: 

f Zx(l-±x 2 ) + 0(x 5 ), x<l 
a(x)*\ y ( 55 ) 
\ l\nAex + 0(^), x>l 

and 

[ fa?(l + ±x 2 ) + 0(x 5 ), x<l 

,H 1f + ^ + *),»i (56) 

where A is defined below Ea. (|35|l . e is the exponential constant and C( x ) ~ Riemann's zeta- function. Combining these 
asymptotics with those given by Eas. 134135(1 . we find 

f S 2^ T > 9 > 9\ 
ml ^{ \ 2 (57) 

I BCS re S ime 

Notice that m z is continuous across T c , however as follows from Eas. (|32l49l5U[l in the AFE phase m z (T,m x ) — 
m z (T,0) + 0(m x ), so its derivative with respect to temperature has a finite discontinuity at T c . 

B. T = 

At zero temperature the parameters t n and rj n are given in terms of the complete elliptic integral of the first and 
second kind, such that 



E(l - A 2 ) - A 2 K (l - A 2 ) 

1 - A 2 
K(l - A 2 ) - E(l - A 2 ) 



*n = " ( 58 ) 

A 

(59) 



1 - A 2 

According to Ea. ((30|l at T = 0: to 2 + m 2 = |, so we can establish the range within which A can lie [cf. Ea.l |29b )] 

< m x < - < A < j < 1 (60) 

where the last strong inequality is just the condition for Hamiltonian's couplings we are implying from physical 
grounds. Note that A max = e/J is the absolute maximum A can reach. Then we can safely make expansions up to 
terms 0(A 2 lnA) 

t n « 1 (61) 

Vn « (62) 
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which can be used in Ea. (|32[l for any g. At large Ising couplings g ^> g\, m z is small, while m x is close to its maximal 
value (A w 2em x /J): 



2A 4e 2 4J 
9 + — In — 

7r 7tJ ee 

1 2 

m„ 

2 



For smaller couplings g < g\ the order parameter to^ is small (A 



5A 



^exp[ - ^(5A - 5)] , BCS regime 



while 



- — to; 
2 



(63) 
(64) 

(65) 
(66) 



is close to its maximal value. The dimerization parameter S to leading order reads 



.9 = 3a 



(67) 



" J 2 ^ exp [ - ^{g x - g)] , BCS regime 



C. Spin susceptibility 

The zero-field spin susceptibility per rung (i.e., per spin) is given by 

d(s z ) _ d 2 f s p fi dip 



Xs 



dh 



M 



dh 2 M 



h M =o 2?r J cosh 2 $£{ip) 



(68) 



Its temperature dependence is similar to that of the dimerized Heisenberg or XY spin chains4i«». In the ordered 
AFE phase it shows the spin-gap behavior. At low temperatures the asymptotics of the integral (|68|l is given by the 
following expression 



1 ^2A SG > 

Xs = 



1 


r Ago] 


J 2 exp 


T \ 



J + Xm 2 . V 7rT 

The characteristic energy scale parameter Age defined as 

A SG = 2em x , (70) 

is natural to call the spin gap. Such definition coincides with the one stated in terms of a magnon4^ The latter is 
equal to the minimal energy needed to flip one spin, and according to Eas. (|19l2U|) . to the one-particle fermion gap. 
The simple relation (|70|) between the spin gap and the alternating part of the exchange b — 2em x [cf. Ea.H26tO] holds 
because the XY model is a free-fermion problem. Accounting for fcrmionic interactions in the XY Z model breaks 
(|70|l already at the mean-field levels Bosonization treatment of the interacting Jordan- Wigner fermionsi^ reveals 
the power-law relationship Ago oc 6 2//3 times some log correction^ 

It is interesting to compare the ratio of the zero-temperature spin gap Ag G and the critical temperature in different 
regimes of couplings. At large g ^> g\ this ratio [cf. Eas. l|63l64l40|l ] 

(71) 

T c g 

is small according to the assumption 10} on the range of the model parameters we are working with. Upon decreasing 
g the ratio increases. It subtly involves interplay of the model couplings, and it is not easy to get its analytic form. 
At the critical coupling g = g* (when A = 0) it can be found from Eas. (|45l65|) . and roughly 



A 



QVlMi ,72) 
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More accurate evaluations of Eas. Q45l65|l . as well as direct numerical calculations show that within the parameter 
range the ratio Ag G /T c does not exceed 1. In the BCS region (g < g\) the ratio is maximal, and moreover, it is 
the universal constant [cf. Eas.l j40l65[l ] 

A° 7T 

— |G = _ ra 1<76j BCS regime ( 73 ) 

1 c 

of any BCS-type theory^! Also, as we can see from Eq.JSJJl, the temperature dependence of the spin gap in this 
regime is exactly as that of the superconductivity gap provided by the BCS theory4£ 

A SG w 3.06T c r5 r < 1, BCS regime . (74) 

D. Charge (pseudospin) susceptibilities 

We will be interested in the pseudospin (i.e., charge) susceptibilities with respect to two external (electric) fields 
which have the same spatial dependence as the Weiss fields h\ nn . Namely, the field along z pseudospin "direction" is 
constant E^ n = — E z , while the field E^ n = (—l) m+n+1 E x is staggered. We define the charge susceptibilities as 



: (i, j) = x or z (75) 

K,=0 



c = dm,, 
Xij ~ dEj 

These quantities can be calculated from equations of Section ITTT1 with the Weiss fields shifted as hi i— > hi + Ei. 

1. T>T C 

After a straightforward but somehow lengthy algebra we obtain for the disordered paraelectric (PE) phase T > T c : 



X% z = XL = (76a) 
/3 1 - 4m 2 z 

4l-M(l-4m2)(t n + Wmj 



XL = T -, W7, : —. \ ~ : ^ ( 76b ) 



(76c) 

l-m z (g+ j(J+ £ Am .) ??n- —t n ) 



where 



dt n (A,(3) 
dp 



Note that even in the limit E z — ► 0, x zz 1S m f &c t the susceptibility under the applied field fl [cf. Eq.QJ] along this 
direction. It does not show critical behavior. On the contrary, from the comparison of Eas. (|32l76b ) we immediately 
see that x% x diverges when T — ► T+. At large g ^> g\ from Eas. (|4QI43|l we find 

Xt z « ^ (78a) 
c _ 4T C 1 

Xxx 7^4 M ( } 

so Xxx diverges with the mean-field value of the critical index. The same divergence can be obtained from the 
analytical treatment in the BCS regime g < g\ for T > T c -C 1 where 

X° zz « ^exp[-^] (79a) 

c 7rJ 7rJ 1 (7QU\ 

X ™ ~ 4e 2 ln£ ~ 4?R 

At these low temperatures one probes the ^-generated charge gap A^g = 1 + A/tt = g\/2 in Xzzj renormalized 
(comparatively to its bare value 1) via the coupling A by the dimer operator average t [cf. Ea. (l24ll ]. 
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2. T <T C 

We will not give here the closed expressions for the four susceptibilities below T c , since there are too cumbersome. 
The important properties of these quantities are the following: xtz ^ s continuous across the transition, and x% x x V 7 " 
is divergent when approaching T c from below. In the regimes of strong Ising coupling and BCS, x%x can be calculated 
analytically, showing that 

X C XX (T -> T~) = \xlAT - T+) (80) 

as one should have expected in the mean-field theory. The transverse pseudospin susceptibilities are "diamagnetic" 
(indicating that the system tends to preserve the psedospin vector's length) and show weaker divergence: 

i ( 9 \ 2 _1 

V 2 (81) 
— ^^JT± T -- BCS regime 

2^/l4C(3) e ' to 

E. Specific heat capacity 

The easiest way to obtain the specific heat capacity in our approach, is to calculate it directly from the mean-field 
energy f^Bf 

(82) 

(here, in fact c = c v ). As usual, at the second-order phase transition temperature the specific heat undergoes a finite 
jump. For two regimes of couplings c can be easily calculated analytically. At high temperatures T>1 above T c 

c«i^(l + Y) (83) 

As one can see from the above equation, spins and pseudospins (m z ) contribute to the specific heat in the same 
fashion. At large g 3> g\, T c is also large. For this regime we find 

q 2 - 4 

9 > <?a (84) 



where c ± = c(T — » T^). Low temperatures (T <C 1) within the disordered PE phase are attainable only at smaller 
couplings g < g\. In the BCS regime (T > T c ) we get 



3J 



c « ^4 + C(T 3 ), BCS regime (85) 



Since pseudospins m z are gapped, at these low temperatures their contribution to the specific heat is exponentially 
suppressed as C(exp[— IS.qq/T\) 1 thus only spins contribute. The latter contribution is equivalent to the one of a 
degenerate free electron gas. Eas. H30|) for the AFE phase in the BCS regime gives 

d -^ = - d -§ + 0(e-^) (86) 

Using it along with asymptotics I|34I35I40I55I56I57|| . we obtain 

4tt T 

°~ " ° + ~ 7C(3) T BCSregime ( 87 ) 

Then the ratio 
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is another universal constant known from the BCS theorycsi Note that the origin of this universality can be traced 
from Eq. I|86|) , which results in that only gapped fermions contribute to the heat capacity in the vicinity of T c in the 
BCS regime. 

Since the specific heat jump is positive at any regime of couplings, it indicates the free energy (|31[) decrease in the 
ordered AFE phase /afe — /pe — m t- 

Upon approaching zero temperature the specific heat capacity decreases exponentially for any coupling g. The 
decay scale is defined by Age which is the smallest gap in the model (the zero temperature charge gap varies from 
Acg = .9a/2 in the BCS regime to Acg = in the regime of large g), so c oc exp[— A$g/T] similar to the spin 
susceptibility. 

VI. SUMMARY AND DISCUSSION 

We study the spin-pseudospin model comprised from the Ising Model in Transverse Field (IMTF) Hamiltonian for 
pseudospins coupled to the spin-^ XK-Hamiltonian, on a triangular lattice. The effective Hamiltonian, similar to 
ours, appears from analyses of the quarter-filled ladder compound NaV205, which can be described by a set of spin 
and pseudospin variables, the latter representing the charge degrees of freedom. Following the proposal of Mostovoy 
and Khomski: 22 - 25 for the scenario of the phase transition in NaA^Os, we include into the Hamiltonian a specific 
inter-ladder spin-pseudospin coupling term allowed by symmetries, linear over the pseudospin operators. 

In the framework of the molecular-field approximation we explore the phase diagram of the model and find that it 
possesses two phases: the disordered phase without charge order or spin gap; and the low-temperature phase containing 
both the anti- ferroelectric (zigzag) charge order and spin gap (without long-range spin order). Such ordered phase is 
experimentally observed below the transition temperature T c = 34 K in NaA^Os. The phase transition in our model 
is of the second kind. We calculate various physical quantities, as the order parameter, spin and pseudospin (charge) 
susceptibilities, specific heat, and find their (mean-field) temperature behavior near T c . 

Our analysis reveals an important property of the phase diagram, not envisaged in the original Mostovoy-Khomskii 
proposal: the inter-ladder spin-pseudospin coupling (e) not only creates the spin dimerization (spin gap) triggered 
by the charge ordering in the region of Ising's coupling where the IMTF can order (<? > gx), but it also generates 
simultaneous appearance of the charge order and spin gap in the case when the IMTF is always disordered (g < g\). 
So, the coupling e destroys the IMTF quantum critical point at g = g\ and results in a continuous evolution of 
T c (g) 7^ into the region g < g\. In that region T c has the exponential BCS-like dependence on model's couplings. 

An interesting feature of model's phase diagram is that near the IMTF critical coupling g\ it shows regimes with 
a re-entrance (for couplings A 7^ 0, e — 0) or with a double re-entrance (for couplings A 7^ 0, e 7^ 0). How pronounced 
the re-entrant behavior is, depends on the interplay of couplings, but the re-entrance is absent when A = 0. Usually, 
re-entrance is due either to disorder effects^ or due to the stabilization of the ordered phase via the release of 
entropy by secondary degrees of freedom^ The IMTF does not show any re-entrance at the mean-field level, 34-36 
but accounting for the first-order fluctuation corrections revealed some re-entrance in the 2D IMTFiSi In that study 
it was assumed to be an artifact of poor approximation. In our case, we have given some qualitative arguments on 
how the re-entrant behavior can be understood from the competition of model's energy scales. However, we cannot 
exclude that the found re-entrance is a mean-field artifact. A more detailed study on this re-entrance will be presented 
elsewhere^ 

Now we proceed with the discussion on the applications of our results to NaA^Ov For NaA^Os we take the 
following datai 17 i 21 i 33 ft = 700 meV, T c = 34 K, A SG = 106 K (A§ G /T C w 3), J H + jJst = 50 meV. The microscopic 
calculations^ give e = 15 meV. Using these data as model's dimensionless parameters [cf. notations 10] e = 0.021, 
J = 0.065, A = 0.026 we find that T c = 0.004 lies on the curve T c (g) in the BCS-regime region with g = 1.99 (J/ = 1.4 
eV), close to the critical coupling g\ — 2.02. The estimates for Jj in NaV205 which is proportional to the in-ladder 
Coulomb repulsion between neighboring rungs, vaiy^i^^ but do not exceed Jj = 1.5 eV££ The whole curve T c (g) 
with the above parameters J, A, e looks very much similar to that shown in Fig. [3for J = 1, A = 0, e = 0.1, without 
appreciable re-entrance. (In any case g = 1.99 < <? m i n .) We would not insist that our mean-field analysis gives the 
quantitatively correct evaluation of T c in terms of the microscopic parameters for NaA^Os, since the uncertainties in 
their values left us with a freedom to make this reasonable fit. However, we think we have grasped the qualitatively 
correct physics. The question why T c is so low in NaV205, or what is the scale that gives such a low T c , can be 
answered as follows: g for NaA^Os lies in the proximity of the quantum critical point of the IMTF (most likely on the 
"disordered side, i.e., g < g\). T c is determined by the two scales: g — g\ and the inter-ladder dimerization coupling 
e. The latter destroys the IMTF quantum critical point and makes the ordering possible at g < g\. 

A flaw of the present study is that we cannot account for the large BCS ratio Ag G /T c observed in NaA^Os. The 
small value of T c indicates that we are either in the BCS regime described by the analytical result l|40|) . or very close 
to it, so the ratio cannot exceed the universal BCS value l|73(l . This gives the spin gap almost two times smaller than 
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the experimental value. This clearly indicates that for more realistic calculations the spin part of the Hamiltonian 
should be modified, and instead of the XY spin model the full three component Heisenberg Hamiltonian should be 
considered. As we have already explained it in Section IV CI this will change results for the BCS ratio. At this point 
we cannot say whether this will be enough in order to explain the experimentally observed magnitude of the ratio. 
Another possibility is of course, the role of modes, neglected in our model, like, e.g., phonons, which can affect the 
BCS ratio. However, we will not speculate more on this point. 

Experiments of various kinds 21 unequivocally demonstrate that what occurs in NaA^Os at T c = 34 K is a (thermal) 
continuous structural phase transition. Also rather large regions of structural (charge-ordering) two-dimensional 
fluctuations are seenS on the both sides of T c . The order parameter critical inde»^i*^iS j3 w 0.17 — 0.19 is close 
to what one expects for a two-dimensional phase transition with a one-component order parameter. So we conclude 
that the universality class of this transition is the 2D Ising. In order to obtain a more complete description of the 
transition in NaV20s, thermal fluctuations in the IMTF-sector of the spin-pseudospin Hamiltonian must be taken 
into account. According to what we said above, we expect the Ising coupling for the NaV2 5 -Hamiltonian to lie near 
the quantum critical point of the 2D IMTF. From the experimental results we expect that physics of NaV2C>5 near 
T c is controlled by the region of the 2D IMTF phase diagram^, where 2D Ising (thermal) fluctuations dominate. 
Although at lower temperatures the change of the nature of these fluctuations into 3D Ising (quantum) may also 
become important. We believe that apart from the applications to NaX^Os, a study of the fluctuations in proximity 
of this quantum critical point is a very interesting problem of its own. 

It is worth pointing out that the pseudospin sector of our Hamiltonian is only quasi two-dimensional, in the sense 
that pseudospins on neighboring ladders are coupled only via the dimcrization coupling e. Accounting for the Ising 
coupling between ladders, as one can see from the triangular lattice mapping shown in Fig^ will cause a frustration. 
The neglect of this coupling for the applications to NaV2 05 is justified from microscopic grounds, but from a general 
point of view it would be interesting to study a Hamiltonian, similar to ours, with a truly two-dimensional IMTF 
part. 
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FIG. 1: Two-dimensional triangular lattice of coupled ladders. A vertical line denotes a single ladder where a dot denotes its 
rung. At each rung (m,n) there is pseudospin T mn and spin S mn (not shown). In the region encirled by the dashed line we 
show four sites involved in each term of the sum in the Hamiltonian Spins and pseudospins residing on two sites of the 
same m-th ladder are coupled via the exchange interaction terms, while two pseudospins from (m — l)-th and (m + l)-th ladders 
are coupled via the dimerization interaction constant e. The pseudospin T£ n ordering pattern shown in the figure corresponds 
to the ordered anti-ferroelectric (zigzag) phase. 
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FIG. 2: Critical temperature of the AFE phase transition as a function of the Ising coupling g at different values of A,£. 
The critical couplings g* — 2 and g\ = 2.6366. Large empty circles on the curves at e — 0.1 indicate the right boundary for 
the BCS-region described by the exponential dependence Q40p. At large values of g (not shown) all curves T c (g) approach the 
asymptotic line T c — g/A. 



